The Newmark Integration Method for Simulation of Multibody Systems: Analytical Considerations

When simulating the behavior of a mechanical system, the time evolution of the generalized coordinates used to represent the configuration of the model is computed as the solution of a combined set of ordinary differential and algebraic equations (DAEs). There are several ways in which the numerical solution of the resulting index 3 DAE problem can be approached. The most well-known and time-honored algorithms are the direct discretization approach, and the state-space reduction approach, respectively. In the latter, the problem is reduced to a minimal set of potentially new generalized coordinates in which the problem assumes the form of a pure second order set of Ordinary Differential Equations (ODE). This approach is very accurate, but computationally intensive, especially when dealing with large mechanical systems that contain flexible parts, stiff components, and contact/impact. The direct discretization approach is less but nevertheless sufficiently accurate yet significantly faster, and it is the approach that is considered in this paper. In the context of direct discretization methods, approaches based on the Backward Differentiation Formulas (BDF) have been the traditional choice for more than 20 years. This paper proposes a new approach in which BDF methods are replaced by the Newmark formulas. Local convergence analysis is carried out for the proposed method, and step-size control, error estimation, and nonlinear system solution related issues are discussed in detail. A series of two simple models are used to validate the method. The global convergence analysis and a computationalefficiency comparison with the most widely used numerical integrator available in the MSC.ADAMS commercial simulation package are forthcoming. The new method has been implemented successfully for industrial strength Dynamic Analysis simulations in the 2005 version of the MSC.ADAMS software and used very effectively for the simulation of systems with more than 15,000 differential-algebraic equations. 1 GENERAL CONSIDERATIONS ABOUT THE NEWMARK METHOD The Newmark method [1] is by far one of the most widely used integration method in the structural dynamics community for the numerical integration of a linear set of second Order Differential Equations (ODE). This problem is obtained at the end of a finite element discretization. Provided the finite element approach is linear, the equations of motion assume the form Mq̈+Cq̇+Kq = F(t) (1) The p× p mass, damping, and stiffness matrices, M, C, and K, respectively, are constant, the force F ∈ Rp depends on time t, and q∈Rp is the set of generalized coordinates used to represent the configuration of the mechanical system. The attractive attributes associated with the Newmark method are: (a) the resulting second order ODEs that govern the time evolution of the system do not have to be reduced to first order, which leads to simpler implementation and a smaller dimension problem; (b) good stability properties and ability to adjust the amount of damping 1 Copyright c © 2005 by ASME introduced into the system; (c) the method has been tested and validated in a vast array of applications spanning many engineering fields. The Newmark family of integration formulas depends on two parameters β and γ: qn+1 = qn +hq̇n + h2 2 [(1−2β) q̈n +2βq̈n+1] (2a) q̇n+1 = q̇n +h [(1− γ) q̈n + γq̈n+1] (2b) These formulas are used to discretize at time tn+1 the equations of motion (1) Mq̈n+1 +Cq̇n+1 +Kqn+1 = Fn+1 (2c) Note that based on Eqs. (2a) and (2b), qn+1 and q̇n+1 are functions of the acceleration q̈n+1, which in Eq. (2c) remains the sole unknown quantity that is eventually computed as the solution of a linear system. This method is implicit and A-stable (stable in the whole left-hand plane) [2] provided [3] γ≥ 1/2 β≥ ( γ+ 1 2 )2